In [1]:
# Bibliotecas que vamos a usar 
import numpy as np
import matplotlib.pyplot as plt
import cv2

Leer la imagen¶

In [2]:
img_path = 'BreCaHAD\BreCaHAD\images\Case_1-01.tif'
groudtruth_path = 'BreCaHAD\BreCaHAD\groundTruth\Case_1-01.json'
In [3]:
# Leer imágenes 
img1 = cv2.imread(img_path)

# Cambiar de BGR a RGB
img1 = cv2.cvtColor(img1, cv2.COLOR_BGR2RGB)
# Mostrar la imagen
plt.imshow(img1)
Out[3]:
<matplotlib.image.AxesImage at 0x1d9119c5150>
In [4]:
#Transformar a escala de grises para luego
img1_gray = cv2.cvtColor(img1, cv2.COLOR_RGB2GRAY)

Extraer canal azul¶

In [5]:
blue_channel = img1[:, :, 0]  # 0 corresponde al canal azul

# Crear una imagen en negro con el mismo tamaño que la original
blue_img = np.zeros_like(img1)

#Asignar el canal azul a la imagen
blue_img[:, :, 0] = blue_channel  
# establecer los canales verde y rojo en 0 (negro)
blue_img[:, :, 1] = 0 
blue_img[:, :, 2] = 0 

blue_img = cv2.cvtColor(blue_img, cv2.COLOR_BGR2GRAY)

plt.imshow(blue_img, cmap = 'gray')
Out[5]:
<matplotlib.image.AxesImage at 0x1d911beb3a0>
In [6]:
# Mínimo y máximo en imagen en escala de grises 
# blue_img.min()
blue_img.max()
Out[6]:
29
In [7]:
img_nuclei = blue_img.copy()
img_nuclei[img_nuclei < blue_img.max()/2 + 1] = 0 # menos a 15 en negro
img_nuclei[img_nuclei >= blue_img.max()/2 + 1] = 255
In [8]:
plt.imshow(img_nuclei, cmap = 'gray')
Out[8]:
<matplotlib.image.AxesImage at 0x1d911c66bf0>
In [9]:
# Para ver a color las partes seleccionadas
img_selected = img1.copy()
img_selected[(img_nuclei == 255) == True] = 255
plt.imshow(img_selected) # Básicamente no hay eosina
Out[9]:
<matplotlib.image.AxesImage at 0x1d911aae2c0>
In [10]:
# Necesitamos que los núcleos sean blancos y el fondo negro
nuclei_foreground = img_nuclei.copy()
nuclei_foreground[nuclei_foreground == 255 ] = 100
nuclei_foreground[nuclei_foreground != 100 ] = 255
nuclei_foreground[nuclei_foreground == 100 ] = 0

"Foreground"¶

In [11]:
# Este es el "sure foreground" que se usa en la transformación de Watershed
plt.imshow(nuclei_foreground, cmap = "gray")
Out[11]:
<matplotlib.image.AxesImage at 0x1d911b25cf0>
In [12]:
# Hay muchos espacios a los interno del núcleo, vamos a tratar de llenarlo con operaciones morfológicas 
# Vamos a usar la operación CLOSING para llenar espacios
kernel = np.ones((6,6),np.uint8)
nuclei_foreground_closing = cv2.morphologyEx(nuclei_foreground, cv2.MORPH_CLOSE, kernel, iterations = 1)
In [13]:
plt.imshow(nuclei_foreground_closing, cmap = "gray")
Out[13]:
<matplotlib.image.AxesImage at 0x1d9129bd7e0>

Background¶

In [14]:
# Vamos a a rellanar "más" los espacios para evitar que haya "background" en los nucleos
kernel = np.ones((3,3),np.uint8)
#background = cv2.morphologyEx(nuclei_foreground_closing, cv2.MORPH_CLOSE, kernel, iterations = 5)
background = cv2.dilate(nuclei_foreground_closing, kernel, iterations = 5)
In [15]:
plt.imshow(background, cmap = "gray")
Out[15]:
<matplotlib.image.AxesImage at 0x1d912a31060>
In [16]:
prueba = img1.copy()
prueba[(background == 0) == True] = 255
plt.imshow(prueba)
Out[16]:
<matplotlib.image.AxesImage at 0x1d912a98700>

Regiones desconocidas¶

In [17]:
unknow_area = cv2.subtract(background, nuclei_foreground)
In [18]:
plt.imshow(unknow_area, cmap = "gray")
Out[18]:
<matplotlib.image.AxesImage at 0x1d91579bd90>

Transformada Watershed¶

In [19]:
# Obtener marcadores
ret, markers = cv2.connectedComponents(nuclei_foreground)
# Agregamos 1 para que ningún marcador sea 0
markers = markers+1
# Marcamos las regiones desconocidas con 0
markers[unknow_area == 255] = 0
In [51]:
nuclei_segmented = img1.copy()
markers = cv2.watershed(nuclei_segmented, markers)
nuclei_segmented[markers == -1] = [255,0,0]
In [52]:
plt.imshow(nuclei_segmented)
nuclei_segmented = cv2.cvtColor(nuclei_segmented, cv2.COLOR_RGB2BGR)
cv2.imwrite('nucleo_water6x6x1.png', nuclei_segmented)
Out[52]:
True

Núcleos segmentados¶

In [22]:
# Obtenemos los núcleos dentro de los bordes 
mask_nuclei = np.zeros_like(img1_gray, dtype=np.uint8)
mask_nuclei[markers != 1] = 255
nuclei_after_water = cv2.bitwise_and(img1, img1, mask=mask_nuclei)
plt.imshow(nuclei_after_water)
Out[22]:
<matplotlib.image.AxesImage at 0x1d91586f340>
In [23]:
nuclei_white = cv2.cvtColor(nuclei_after_water, cv2.COLOR_RGB2GRAY)
# fondo negro y objetos en blanco
nuclei_white[nuclei_white != 0] = 255
plt.imshow(nuclei_white, cmap = "gray")
Out[23]:
<matplotlib.image.AxesImage at 0x1d9156ae830>

Componentes conectados¶

In [24]:
# Parece que el algoritmo que esta usando connectedComponentsWithStats() considera que los objetos en bordes opuestos
# son los mismos. Para correguirlo, vamos a eliminar una pequeña sección de los bordes 


# Eliminar objetos pequeños en los bordes
border_size = 5  # Tamaño del borde para considerar
nuclei_white[:border_size, :] = 0
nuclei_white[-border_size:, :] = 0
nuclei_white[:, :border_size] = 0
nuclei_white[:, -border_size:] = 0
In [25]:
num_labels, labels, stats, centroids = cv2.connectedComponentsWithStats(nuclei_white, connectivity = 8)
In [26]:
# Cantidad de componentes conectados
num_labels
Out[26]:
383
In [27]:
# Vemos el área de los objetos
stats[:, cv2.CC_STAT_AREA]
Out[27]:
array([975024,    169,  13316,    211,     22,    254,    181,      4,
         4459,    159,     13,   1212,    699,    686,    863,     20,
           91,     64,    205,     22,    559,     87,    176,    115,
           19,     59,    734,   1273,    171,    487,    355,    405,
          575,    765,    107,    667,     30,    655,    157,    430,
          632,    498,    837,    480,    644,    268,   2704,    246,
          755,     45,    483,    508,    118,   1598,    483,    462,
         1122,    602,    560,    335,   1288,    254,     35,    136,
         4333,   1177,   4355,    358,     66,   1257,   4493,     14,
           36,    442,   1781,   2106,   9330,     22,    238,    998,
         1469,    227,    343,     95,   1267,   7620,   1047,    218,
           20,    373,     64,     93,   6743,     34,    476,     46,
         1579,   1359,     40,    137,    721,     85,    706,     41,
            3,    353,    774,    951,    351,    112,  11965,     47,
          746,    300,     31,   5400,    441,    728,   1142,    383,
          102,    122,    350,   1568,     81,   3471,    447,      5,
          138,   5452,   4353,    266,    501,    731,    139,    439,
         1613,     39,     40,    573,     14,     39,   1115,     48,
          830,    142,     71,    583,    843,     31,    106,     55,
          243,    327,   8312,    279,   1025,    270,     38,    246,
          658,     61,    654,     40,   4546,    101,    299,    343,
           10,    935,    401,    323,     25,    216,   5556,    102,
          303,    435,   1180,     22,     28,    362,   1360,     62,
         2346,    350,     32,   5128,    264,    176,    150,   6282,
          126,    950,   1691,    533,     20,   1146,    444,    560,
           56,   1273,   2659,    457,   2376,     60,     11,   2352,
           21,     50,   2165,    448,   2073,    140,   4947,   1794,
         6120,    112,     91,     23,     86,    237,   1920,    245,
          587,   6102,     47,   4806,   1072,     94,     37,     12,
           89,   1100,    668,     19,    553,    533,    607,    114,
         8263,    169,    558,   2546,   4412,   2204,     66,     34,
          993,    184,    582,   2076,    195,    111,   1522,    992,
         2262,   1920,   1055,    843,    753,    518,   8376,    679,
         2315,     21,     54,    682,    889,   1772,    672,     14,
          864,   1415,   9791,   2640,   1343,   2534,     37,    794,
          835,    197,     14,     10,    969,   1223,    641,   1152,
          520,     27,     18,    445,   1122,   1594,     81,   1007,
          339,   1471,     34,   8403,    454,   1179,    160,    159,
          496,   1277,  11703,     24,    773,   1813,   2974,    293,
           85,     69,    150,   2908,    724,    149,   5354,   1295,
          527,   1131,     36,    983,    743,     44,     47,   5025,
         2956,    107,    162,    572,     32,   1620,    995,    980,
          462,    227,    610,   1171,    155,    181,   3013,   3007,
         1972,   1421,    548,    266,     92,   5830,    134,   1023,
          306,     90,    708,     16,    887,    275,   1741,     21,
          290,    678,    134,    767,    857,   1279,    635,    183,
         1409,   1434,    381,     40,    370,     86,    674,   1030,
          659,    416,    351,    739,   1043,     92,    134],
      dtype=int32)
In [28]:
# Tambien se pueden extraer otras caracteristicas de esos componentes, como:
# Altura del cuadro delimitador 
# stats[:, cv2.CC_STAT_HEIGHT] 

# Ancho del cuadro delimitador
#stats[:, cv2.CC_STAT_WIDTH]

# Los centroides están en variable `centroids`

Retirar componentes conectados pequeños¶

In [29]:
# Definimos un tamaño mínimo
min_size = 800 # ELiminamos los objemos menores a este tamaño

# Generamos un "lienzo" en negro
nuclei_without_smalls = np.zeros((nuclei_white.shape))

for i in range(1, num_labels):
    if stats[i, cv2.CC_STAT_AREA] >= min_size:
        nuclei_without_smalls[labels == i] = 255

plt.imshow(nuclei_without_smalls, cmap = 'gray')
Out[29]:
<matplotlib.image.AxesImage at 0x1d9158013c0>
In [30]:
nuclei_with_borders = nuclei_without_smalls.copy()

# Ponemos en negro los bordes obtenidos con Watershed
nuclei_with_borders[markers == -1] = 0
plt.imshow(nuclei_with_borders, cmap = 'gray') # Lo que se gana es poco
cv2.imwrite('nuclei_gray.png', nuclei_with_borders)
Out[30]:
True

Resultado¶

In [31]:
# Núcleos 
nuclei_color = img1.copy()
nuclei_color[(nuclei_without_smalls == 0) == True] = 0
plt.imshow(nuclei_color)
nuclei_color = cv2.cvtColor(nuclei_color, cv2.COLOR_RGB2BGR)
cv2.imwrite('nuclei_color.png', nuclei_color)
Out[31]:
True
In [32]:
# Fondo
background_color = cv2.subtract(img1, nuclei_color)
plt.imshow(background_color)
Out[32]:
<matplotlib.image.AxesImage at 0x1d915ff7e50>

Métricas¶

Leer groundTruth¶

In [33]:
import json

with open(groudtruth_path, 'r') as file:
    # Carga el contenido del archivo JSON
    groundTruth_img1 = json.load(file)
In [34]:
 # Eliminar la clase lumen y non_lumen que forman parte del background
del groundTruth_img1['lumen']
del groundTruth_img1['non_lumen']
In [35]:
img1.shape
Out[35]:
(1024, 1360, 3)
In [36]:
img1_points = img1.copy()
nuclei_points = nuclei_color.copy() 
height = img1.shape[0]
width = img1.shape[1]
width
Out[36]:
1360
In [37]:
# Asignar un color para cada clase 
cell_color = {
    'mitosis': (0, 0, 255),        # azul
    'non_tumor': (0, 255, 0),    # Verde
    'apoptosis': (255, 0, 0),      # rojo
    'tumor': (255, 255, 0),         # Amarillo
    'non_mitosis': (0, 0, 0),    # Negro
}

Dibujar y determinar cantidad de núcleos identificados¶

In [38]:
# Dibujar círculos en la imagen basándote en las coordenadas del JSON
img_with_points =  np.zeros_like(img1)

# extraer valores
total_nuclei_groundTruth = 0
identified_nuclei = 0

for cell_type, nuclei in groundTruth_img1.items():
    color = cell_color.get(cell_type, (255, 255, 255))  # blanco si no está definido
    for point in nuclei:
        x = int(point['x'] * width)
        y = int(point['y'] * height)  # Usar img.shape[0] para obtener la altura
        cv2.circle(img1_points, (x, y), 10, color, -1)  # Dibuja círculos de radio 10 en color rojo (BGR)
        cv2.circle(nuclei_points, (x, y), 10, color, -1)
        cv2.circle(img_with_points, (x,y), 10, (255,255,255) , -1)
        # Cantidad total de núcleos
        total_nuclei_groundTruth +=1
        # Cantidad de ocaciones donde el pixel del "groudtruth" corresponde a un pixel de núcleo
        pixel = nuclei_color[y][x]
        
        if (pixel[0] != 0 and pixel[1] != 0 and pixel[2] != 0):
            identified_nuclei +=1
            
print("La cantidad total de núcleos en el `groundTruth` es", total_nuclei_groundTruth)
print("De eso se segmentaron ", identified_nuclei)
print ("Para un procentaje de", identified_nuclei/total_nuclei_groundTruth*100)
# Este porcentaje es igual a la sensibilidad: 
    # verdaderos_positivos / (verdaderos_positivos + falsos_negativos)
    # identified_nuclei / (identified_nuclei + (total_nuclei_groundTruth - identified_nuclei))
    
        

La cantidad total de núcleos en el `groundTruth` es 141
De eso se segmentaron  128
Para un procentaje de 90.78014184397163
In [39]:
nuclei_points.shape
Out[39]:
(1024, 1360, 3)
In [40]:
# Mostrar resultado 
plt.imshow(img1_points) # los puntos corresponden al groudtruth
Out[40]:
<matplotlib.image.AxesImage at 0x1d91b2a8250>
In [41]:
plt.imshow(nuclei_points)
Out[41]:
<matplotlib.image.AxesImage at 0x1d91b2ef700>
In [42]:
plt.imshow(img_with_points)
Out[42]:
<matplotlib.image.AxesImage at 0x1d91b352cb0>

Regiones segmentadas y no presentes en el "groundTruth"¶

In [43]:
def groundTruth_inside_square(groundTruth, x, y, width, height, image_width, image_height):
    for cell_type, nuclei in groundTruth.items():
        for point in nuclei:
            pixel_x = int(point['x']*image_width)
            pixel_y = int(point['y']*image_height)
            if (x <= pixel_x <= x + width) and (y <= pixel_y <= y + height):
                return True
                break
    return False
In [44]:
nuclei_without_smalls_8u = nuclei_without_smalls.astype(np.uint8) #necesita que sean enteros sin signo
num_labels, labels, stats, centroids = cv2.connectedComponentsWithStats(nuclei_without_smalls_8u, connectivity = 8)

n_nuclei_segmented = 0
not_in_groud_truth = 0

squere_components = np.zeros_like(img1)
# Dibujar el cuadro delimitador de cada componente conectado
for i in range(1, num_labels):  # Empezar desde 1 para omitir el fondo (etiqueta 0)
    x, y, w, h, area = stats[i]
    n_nuclei_segmented += 1
    if (not groundTruth_inside_square(groundTruth_img1, x, y, w, h, width, height)):
        not_in_groud_truth += 1
        
    # Dibujar un cuadrado delimitador en el lienzo
    cv2.rectangle(squere_components, (x, y), (x+w, y+h), (255, 255, 255), -1)
plt.imshow(squere_components)

print("La cantidad de objetos segmentados fue", n_nuclei_segmented)
print("De esos", not_in_groud_truth , "no estaban en el `groudtruth` ")

La cantidad de objetos segmentados fue 126
De esos 42 no estaban en el `groudtruth`
In [45]:
plt.imshow(img_with_points)
Out[45]:
<matplotlib.image.AxesImage at 0x1d91ddcd900>

Precisión / Especificidad¶

In [46]:
# verdaderos_positivos / (verdaderos_positivos + falsos_positivos)
identified_nuclei / (identified_nuclei + not_in_groud_truth)
Out[46]:
0.7529411764705882

Sensibilidad¶

In [47]:
identified_nuclei / (identified_nuclei + (total_nuclei_groundTruth - identified_nuclei))
Out[47]:
0.9078014184397163

Anexos¶

In [48]:
# Métricas para imagenes de prueba fueron
# Case_1-01.tif | Precisión: 0.7555555555555555 |  sensibilidad: 0.9645390070921985
# Case_1-02.tif | Precisión: 0.7978723404255319  |  sensibilidad: 0.8771929824561403
# Case_1-03.tif | Precisión: 0.770949720670391 |  sensibilidad: 0.9019607843137255
# Case_1-04.tif | Precisión: 0.7536231884057971 |  sensibilidad:0.9043478260869565
# Case_1-05.tif | Precisión: 0.7441860465116279 |  sensibilidad: 0.9552238805970149
# Case_1-06.tif | Precisión: 0.7479338842975206 |  sensibilidad: 0.8916256157635468
# Case_1-07.tif | Precisión: 0.871244635193133 |  sensibilidad: 0.8982300884955752
# Case_1-08.tif | Precisión: 0.76 |  sensibilidad: 0.957983193277311